%% estimaStates 

%% 1. Kalman filter and State decomposition 
[stt,etamat,smooth_st,yfor,vstar,logL,countmat,countobs,countszero]=...
    feval(addsol.funcKfilter,parmode_top,funcmod,Y,trainvec,solveopt,addsol); 

%% 3. State decomposition for observables or time aggregated shocks
% Recall POS_DECOMPSTATES and DECOMPSTATES_NAMES are different from ZNAMES
% and POS_OBSINSTATES in that they adjust to plot the cummulators at the
% end
close all;
figvec=zeros(length(znames),1);
% MIXFREQ_VEC indiscator ==1  if vector is mixed frequency 
mixfreq_vec=figvec; 
for ii=1:length(znames);
    if strcmp(znames{ii},decompstates_names{ii})==1
        figvec(ii)=plot_decompositions(Y(:,ii),squeeze(countmat(:,pos_decompstates(ii),:)),znames(ii),shonames,sample,[decompstates_names{ii},': decomposition per shock'],0.2,0);
    else
        if flag_mixfreq==0
            error('Strings must match')
        end
        mixfreq_vec(ii)=1; 
        figvec(ii)=plot_decompositions(Y(posinsamp_lf,ii),squeeze(countmat(posinsamp_lf,pos_decompstates(ii),:)),znames(ii),shonames,indecomp.sample_lf,[decompstates_names{ii},': decomposition per shock'],0.2,0);
    end
end
indecomp.mixfreq_vec=mixfreq_vec; 
cd(outpath); 
for ii=1:length(znames);
    if ii == 1
        print(figvec(ii),'-dpsc','decomposition shock'); 
    else
        print(figvec(ii),'-dpsc','decomposition shock','-append'); 
    end
end
cd(cucd); 
%% Plot per shock 
%pospershock=cellposition({'unemployment','finding','separation','vacancies','gdp','I','C'},stanames); 
pospershock=cellposition(znames,stanames); 

%% Decomposition with Neutral and MEI shocks only 
% [figvec_zphi,decmat_zphi,vrat_zphi]=plot_decompositions_pergroup_mixfreq(smooth_st(:,pos_decompstates),countmat(:,pos_decompstates,:),decompstates_names,shonames,{'Neutral technology','MEI'},'Neutral and MEI shocks only',0.1,0,indecomp);
% cd(outpath);
% print(figvec_zphi,'-dpdf','Neutral MEI only');
% cd(cucd);




%% Decomposition with all but measurement errors

% [figvec_nome,decmat_nome,vrat_nome]=plot_decompositions_pergroup_mixfreq(smooth_st(:,pos_decompstates),countmat(:,pos_decompstates,:),decompstates_names,shonames,{'Neutral technology','MEI','Job destruction','Discount factor','Aggregate demand','Matching'},'All but measurement errors',0.1,0,indecomp);
% 
% indecomp.decompmat_zphi=decmat_zphi; 
% indecomp.vratio_zphi =vrat_zphi; 
% indecomp.decompnames=decompstates_names; 
% indecomp.label_zphi={'Neutral & MEI Only'}; 
% indecomp.actual    =smooth_st(:,pos_decompstates); 
% indecomp.vratio_nome=vrat_nome; 
% 
% cd(outpath);
% print(figvec_zphi,'-dpdf','All but Measurement errors');
% save counterfactuals indecomp smooth_st pos_decompstates countmat; 
% cd(cucd);


figvec=plot_decompositions_pershock(smooth_st(:,pospershock),countmat(:,pospershock,:),stanames(pospershock),shonames,sample,[znames{ii},': decomposition per shock'],0.2,0);
cd(outpath); 
for ii=1:length(znames);
    if ii == 1;
        print(figvec(ii),'-dpsc','decomposition series'); 
    else
        print(figvec(ii),'-dpsc','decomposition series','-append'); 
    end    
end
cd(cucd); 

% Create table with the states and innovations 
cd(outpath); 
save statemat smooth_st stt stanames shonames etamat; 
if flag_mixfreq==1
    outputcell = crtcell(smooth_st(:,1:end-out.nzLF),num2cprec(sample),stanames,{'SmoothStates'}); 
else
    outputcell = crtcell(smooth_st,num2cprec(sample),stanames,{'Smooth States'}); 
end

if ~isunix
    xlswrite('tab states',outputcell,'SmoothSt');
else
    savecell(outputcell,[],outpath,['tab states', '_', 'SmoothSt']);
end
if ~isunix
    xlswrite('tab states',crtcell(etamat,num2cprec(sample),shonames,{'Smooth Innovation'}),'SmoothEta'); 
else
    savecell(crtcell(etamat,num2cprec(sample),shonames,{'Smooth Innovation'}), ...
        [], outpath, ['tab states', '_','SmoothEta']); 
end
cd(cucd); 

%% Plot stochastic shocks and their innovations 
% Shocks_sel variable must be defined in the Settings.m file 
if exist('shocks_sel','var')==1 && ~isempty(shocks_sel)
    stshocks_pos=cellposition(shocks_sel,stanames);
    figvec=zeros( length(shocks_sel) );
    for ii=1:length(shocks_sel);
        figvec(ii)=figure;
        subplot(2,1,1);
        plot(sample,smooth_st(:,stshocks_pos(ii)),'LineWidth',2 );
        title([shocks_sel{ii},': exogenous process'],'FontSize',14);
        xlim([sample(1) sample(end)+0.1]);
        nberdates_ajd('y');
        subplot(2,1,2);
        plot(sample,etamat(:,ii),'LineWidth',2 );
        title([shonames{ii},': innovations'],'FontSize',14);
        xlim([sample(1) sample(end)+0.1]);
        subtitulo(extrsubf(outpath,cucd));
    end
    cd(outpath);
    for ii=1:length(shocks_sel);
        if ii == 1;
            print(figvec(ii),'-dpsc','Shocks');
        else
            print(figvec(ii),'-dpsc','Shocks','-append');
        end
    end
    if ~isunix
        ps2pdf('psfile', 'Shocks.ps', 'pdffile', 'Shocks.pdf', 'deletepsfile',1);
    end
end
cd(cucd); 